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The  principal  objective  of  this  research  was  to  develop  a  method 
for  determination  of  the  response  to  random  excitation  of  structures 
having  geometric  nonlinearities .   This  would  account  for  large  or  finite 
deflections  of  such  structures  whereas  work,  done  in  the  past  has  accounted 
only  for  small  deformations  leading  to  linearized  equations. 

The  method  of  attack  was  to  expand  the  transverse  deflection  of  the 
structure  in  terms  of  the  eigenfunctions  of  the  linear  problem  and  the 
longitudinal  deflections  in  terms  of  orthogonal  functions  satisfying  the 
boundary  conditions.   The  equations  governing  the  series  coefficients  of 
the  expansions  are  then  easily  derived  from  the  Euler-Lagrange  variational 
equations  and  are  shown  to  be  nonlinearly  coupled.   The  assumption  of 
uncorrelated  loading  then  permitted  the  identification  of  the  phase  space 
of  the  series  coefficients  with  a  Markoff  process  and  thereby  permitted 
derivation  of  the  Fokker-Planck  equation  governing  the  first  order 
probability  density  function  of  the  series  coefficients.   A  stationary 
solution  of  this  equation  was  obtained.   While  the  solution  found  is 
valid  for  several  structural  elements  with  various  boundary  conditions, 


VI 


a  simply  supported  beam  and  a  simply  supported  plate  were  investigated 
in  detail.   It  was  found  that  the  type  of  nonlinearities  considered  had 
the  effect  of  reducing  the  mean-squared  response.   Furthermore,  it  was 
found  that  the  mean-sqdared  value  of  the  first  mode  still  represented  a 
good  estimate  of  the  total  mean-squared  deflection  of  the  structure  but 
in  contradistinction  to  the  linear  theory  the  effects  of  the  higher 
modes  must  be  considered  in  calculating  the  response  of  the  first  mode. 
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CHAPTER  I 


INTRODUCTION 


The  rocket  and  jet  propulsion  systems  of  modern  air  and  space 
vehicles  have  given  rise  to  a  host  of  new  problems  in  the  theory  of 
vibrations.   The  pressure  fields  generated  from  these  systems  vary 
randomly  in  both  time  and  space  over  a  wide  range  of  frequencies. 
These  random  pressure  fields  can  cause  severe  vibrations  to  the  vehicles 
and  their  components.   The  importance  of  these  types  of  vibrations  is 
evidenced  by  the  fact  that  in  1958  and  again  in  1963  a  special  summer 
program  on  random  vibrations  was  held  at  the  Massachusetts  Institute  of 
Technology.   The  two  volumes  of  the  book  "Random  Vibrations,"  edited  by 
Stephen  H.  Crandall  (1,  2)  ,  were  outgrowths  of  these  programs.   In  the 
first  volume  it  was  pointed  out  that  "Some  of  the  strictly  mechanical 
problems  .  .  .  are  still  incompletely  understood  and  the  tools  for 
handling  them  are  relatively  crude."   This  statement  stands  in  evidence 
as  to  the  need  for  research  in  this  important  field. 

The  first  satisfactory  treatment  of  stochastic  motion  was  presented 
by  Einstein  who  studied  the  random  motion  of  a  free  particle.   Smoluchowski 
generalized  this  theory  to  other  types  of  Brownian  motion  and  since  then 
many  important  contributions  have  been  made  to  the  theory,  notably  by 
Fokker,  Planck,  Ornstein,  Uhlenbeck,  Chandraseklar ,  Kramers,  and  others. 


Numbers  in  parentheses  refer  to  the  List  of  References, 
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On  the  purely  mathematical  side,  some  of  the  outstanding  contributions 
have  been  made  by  Wiener,  Kolmogoroff ,  Feller,  and  Doob.    Generally 
speaking,  the  theory  of  linear  lumped  parameter  systems,  i.e.,  systems 
governed  by  ordinary  differential  equations,  excited  by  stochastic 
driving  functions  has  reached  a  relatively  high  level  of  sophistication 
and  the  technique  of  spectral  analysis  developed  by  Rice  (5)  is  generally 
adequate  for  the  solution  of  problems  of  this  type. 

For  nonlinear  lumped  parameter  systems,  the  theory  is  not  quite 
as  refined.   For  small  nonlinearities  an  approximation  procedure 
known  as  "equivalent  linearization"  has  been  developed  (6,  7,  8).   The 
procedure  consists  of  replacing  the  nonlinear  system  with  an  equivalent 
linear  system.   Tl^e  crux  of  the  procedure  is  to  choose  the  stiffness 
matrix  of  the  linear  system  so  that  the  mean  square  error  of  the  governing 
equations  is  a  minimum.   When  the  nonlinearities  are  not  small  this 
method  obviously  fails. 

Another  approximate  procedure  applicable  to  systems  with  small 
nonlinearities  is  the  "perturbation  method"  (see  reference  2).   The 
idea  in  this  method  is  to  assume  a  series  solution  in  powers  of  a 
parameter  which  represents  the  nonlinear  component  of  the  system.   A 
set  of  linear  differential  equations  governing  the  coefficients  of  the 
expansion  can  then  be  generated.   These  linear  equations  can  then  be 
handled  by  more  established  techniques.   For  systems  with  nonlinearities 
which  are  not  small  an  alternate  approach  must  be  made. 


For  a  complete  set  of  references  see  entries  (3)  and  (4)  in  the 
List  of  References. 


Such  an  alternate  approach  exists  in  the  identification  of  the 
trajectory  in  the  phase  space  with  a  Markoff  process.   This  assumption, 
together  with  that  of  Gaussian  input  permits  the  derivation  of  a  partial 
differential  equation,  known  as  the  Fokker-Planck  equation,  governing 
the  probability  density  function  of  the  response  of  the  system.   Tliis 
concept  will  be  more  fully  exploited  in  Chapters  IV  and  V.   We  mention 
here  only  its  shortcomings.   These  are  that  the  assumption  of  a  Markoff 
process  implies  white  noise  input  which  can  be  questioned  with  regard 
to  realizability  and  furthermore,  solutions  to  the  Fokker-Planck 
equations  are  known  only  in  a  few  special  cases  (9,  10,  11).   Despite 
these  shortcomings  it  appears  to  be  the  most  fruitful  approach  for 
analyzing  nonlinear  systems  excited  by  stochastic  functions. 

The  picture  is  worse  for  continuous  systems,  i.e.,  systems  governed 
by  partial  differential  equations.   The  linear  theory  has  been  attacked 
by  Eringen  (12)  who  used  the  generalized  Fourier  analysis  developed  by 
Wiener  (13)  to  the  problems  of  vibrating  beams  and  plates.   Other  special 
linear  continuous  systems  have  been  studied  by  many  authors  (14-20). 
Chapter  III  contains  an  outline  of  the  standard  approach  to  these 
problems. 

As  in  lumped  parameter  systems,  the  method  of  ''equivalent 
linearization"  has  recently  been  used  to  study  continuous  systems 
governed  by  equations  containing  small  nonlinearities  (21).   Until  now 
there  has  been  no  analysis  presented  for  the  problem  of  continuous 
systems  which  are  governed  by  equations  containing  large  nonlinearities. 
It  is  the  purpose  of  this  paper  to  show  how,  with  the  assumption  of  a 
Markoff  process  and  Gaussian  input,  the  Fokker-Planck  equation  can  be 


used  to  determine  the  response  of  nonlinear  elastic  plates  excited  by 
uncorrelated  stochastic  loadings. 

The  nonlinearities  of  the  plate  which  are  considered  are  those 
arising  from  geometrical  considerations.  We  still  make  the  usual 
assumptions  of  small  strains,  shears,  and  extensions  compared  to  unity 
so  that  Hooke's  law  is  valid,  but  we  allow  the  rotations  to  be  moderately 
large,  i.e.,  small  compared  to  unity  but  large  compared  to  the  shears 
and  extensions.   This  will  be  the  result  when  the  deflections  of  a 
plate  are  not  small  relative  to  its  thickness,  a  situation  very  common 
in  structures  subjected  to  vibrations  produced  from  jet  and  rocket 
engines. 


CHAPTER  II 


THEORY  OF  PLATES 


The  theory  of  the  large  forced  vibrations  of  elastic  plates  is 
presented  in  this  chapter.   The  equations  governing  the  large 
vibrations  of  beams  are  obtained  by  properly  simplifying  the  expression 
for  the  kinetic  potential. 

2.1.   Analysis  of  Deformation 

We  start  our  study  of  plates  with  the  geometry  of  deformation. 
Let  the  position  of  a  point  in  an  undeforraed  body  at  time  t^  be 
designated  by  x,  y,  z  referred  to  some  rectangular  coordinate  system 
X,  Y,  Z.   After  deformation,  at  time  t,  the  point  has  moved  to  a  new 
location  "^  ,  7]      >  ^    referred  to  the  Same  X,  Y,  Z  coordinate  system. 
Then  we  can  write 


(2.1) 
The  functions  yCL  ,  V  ,  !>/"  represent  the  projections  onto  the  X,  Y,  Z 
axes  of  the  displacement  at  time  t  of  a  point  which  was  at  the  position 
X,  y,  z  at  time  t^ . 
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If  in  equations  (2.1)  we  set  x  =  Xq  =  constant,  y  =  Yo  =  constant, 
we  obtain  the  equations 


(2.2) 


These  are  the  equations  of  a  line  which  at  time  t   was  parallel  to  the  Z 
axis.   At  time  t  it  is  some  curvilinear  line  as  indicated  in  Figure  1. 
Thus,  while  the  x,  y,  z  form  a  rectangular  coordinate  system  in  the 
undeformed  body,  they  form  a  curvilinear  coordinate  system  in  the  deformed 
body.   Therefore,  if  we  speak  of  a  stress  (J^,   we  mean  the  stress  on  a 


X' 


veto')  ^, 


MCi) 


■Y 


pa) 


Fig.  1.-  Deformation  of 
a  coordinate  line. 


^ 


^^ 


tfzz 


Y 


^% 


Undeformed     Deformed 

Fig.  2.-  Deformation  of  an 
element  of  volume. 

a.  Lagrangian  coordinate  system. 

b.  Eulerian  coordinate  system. 
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surface  which  was  originally  perpendicular  to  the  Z  axis  and  acting  in 
a  direction  normal  to  this  surface  (see  Figure  2).   On  the  other  hand  a 
stress   (iLi,   will  be  the  stress  on  a  surface  which  in  the  deformed  body 
is  perpendicular  to  the  Z  axis  and  acting  in  the  Z  direction  (see  Figure  2) 

Having  dispensed  with  these  preliminary  considerations  of  our 
coordinate  system,  we  may  now  proceed  to  our  analysis  of  strain.   If  we 
introduce  the  following  parameters 


e..  -  ^^- 


>-T^     a  X 


3^       a^ 


^'■'-  '  at 
•^  3      az       5^ 


i^.  =  i-(#-i¥) 


(2.3) 


-8- 


(2.4) 


then  the  components  of  strain  can  be  defined  as  (22) 

£x.  ^  e^M[el  ^iic^^^uS  ^c-t^xz-uy/] 


We  now  make  the  assumption  that  the  strains  and  rotations  of  a 
volume  element  are  small  compared  to  unity.   The  assumption  of  small 
strains  permits  us  to  use  Hooke ' s  Law  which  is  generally  not  valid  for 
large  strains.   When  we  add  the  assumption  of  small  rotations, 
equations  (2.4)  can  be  written  as  (22) 
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(2.5) 


^xx"  ^^y^^T^-^'^^z) 


^^^"  ^3aHCo^^  +  txy/j 


^z.  -^zz^i((x/>tx/p 


^xM  ~  ^"^-y    ijJx^tt 


A  further  consequence  of  our  assumption  is  that  the  quantities  (x/;!  , 
L)Ju  )  CuC  can  be  interpreted  as  the  components  of  the  rotation  vector  of 
a  volume  element.   The  retention  of  their  squares  in  equations  (2.5)  is 
due  to  the  fact  that  while  these  quantities  may  be  small  compared  to 
unity  they  may  be  large  compared  to  the  strains  so  that  their  squares 
may  be  of  the  same  order  of  magnitude  as  the  strains.   For  flexible 
bodies  this  is  quite  often  the  case. 

We  are  now  ready  to  start  our  analysis  of  thin  plates.   We  consider 
that  the  plate  is  of  constant  thickness  h  and  that  the  XY  plane  forms 


For  an  excellent  treatment  of  these  points  see  reference  (23), 
p.  47. 
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the  middle  surface  of  the  plate  before  deformation.   Thus,  the  Z  axis 
is  normal  to  the  plate.   Therefore,  after  deformation  the  xy  surface 
will  form  the  middle  surface  and  the  z  axis  will  be  normal  to  that 
surface . 

As  is  always  the  case  in  a  "Strength  of  Materials"  type  of  analysis 
we  must  assume  a  form  for  the  displacement  field.  For  a  thin  plate  with 
moderately  large  rotations  we  take, 


lLiX,l^,Z,t)   ~    ILL'^,^,±)  -Z. 


3X- 


ucx,i^,z,t)  =  Arc%,u,t) -z  1^  ^^'^^ 


"Ur  L-)C,a  -L^t)  -  AJJ-C-/.,u  ±:) 


As  was  pointed  out  by  Biot  (24)  equations  (2.6)  are  those  of  the 
von  Karraan  plate  theory.  They  are  tantamount  to  assuming  that  a  straight 
line  originally  normal  to  the  middle  surface  of  the  plate  remains 
straight  and  normal  to  the  middle  surface  after  deformation.  This  is  the 
assumption  of  classical  plate  theory.  Furthermore  equations  (2.6)  permit 
a  stretching  of  the  middle  surface  of  the  plate  which  in  classical  theory 
is  omitted. 

The  use  of  equations  (2.6)  yields  the  following  results  for  the 
parameters  given  by  equations  (2.3) 


■11- 


^zz^O 


e^,  =  0 


Ccz  =  0 


^      a 


(2.7) 


The  rotations  in  the  plane  of  the  plate  can  be  considered  fairly  small 
so  that  in  using  equations  (2.7)  to  calculate  the  strain  components  we 


neg 


lect  [x)-,  . 


z 


Tlius,  making  use  of  (2.7)  the  strain  components  become 


■12- 


C.-W.TC 


(2.8) 


£„-0 


=  0 


where 


•X-Tt 


(2.9) 


^^^       3^       a-^      a^  a^ 


The  quantities  ^x-^  ,  €,yu  ,  ^»a   ^^^  obviously  the  components  of  the 
strain  of  the  middle  surface. 
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That  C-xi  =  Cw£  =  0  points  up  the  fact  that  we  are  neglecting 
the  deflections  due  to  the  transverse  shear  stresses,  an  approximation 
that  is  valid  in  most  cases. 


2.2.   Equations  of  Motion 

To  derive  the  equations  of  motion  governing  the  forced  vibrations 
of  plates  we  employ  Hamilton's  principle  in  the  form 


(2.10) 


where  K  is  the  kinetic  potential  of  the  system  given  by 


K  =  T-V 


(2.11) 


Here  T  is  the  kinetic  energy  of  the  system  and  V  the  potential. 
Equation  (2.10)  states  that  the  displacement  field  assumed  by  a  system 
is  such  that  the  kinetic  potential  is  an  extremum. 

For  an  elastic  structure  V  consists  of  two  parts,  the  internal 
strain  energy  Vg  and  the  potential  energy  of  the  external  loads  V  . 
Thus, 


Vs- 


<^-t^Tux  +  *'^4^Hii  "*"  ^z-^zz 


"*"  <^^*^tw+<ry.i€.^i  +■  <^'^z^^i 


(2.12) 


dv 


■14- 


and 


v.- 


^xi+  {^v+K^ur 


cLs 


(2.13) 


where  the  first  expression  is  integrated  over  the  volume  of  the  body 
and  the  second  over  the  entire  surface  of  the  body.   The  quantities 

p  ,  l|_  ,  Pj-   are  the  components  of  the  boundary  traction  in  the  X,  Y, 
Z  directions,  respectively.   Since  in  our  analysis  we  have  assumed  small 
strains,  no  distinction  has  been  made  between  integration  over  the 
deformed  and  undeformed  body  nor  between  the  stresses  and  the  pseudo- 
stresses  a^^^,   which  actually  should  be  used  in  equation  (2.12). 

For  a  plate  with  no  shear  stresses  on  top  or  bottom,  equation  (2.13) 
reduces  to 


V.  -  - 


or 


2  =  X 


IXlds 


■rhii 


\<^  ^|.i^) 


dtdl 


(2.14) 


V.= 


a  (%,\j,il:)iviiL,[j^i]dS 


h/2. 


(I^XL  ^L-^  ^[^uj)diM 


V 


-h/i 


(2.15) 


■15- 


where 


f '^.^.i'  =[f. 


-K/iL 


(2.16) 


Z---K/i 


and  the  line  integral  is  taken  around  the  entire  cylindrical  boundary 
of  the  plate.   We  consider  only  the  case  when  the  forces  on  the 
cylindrical  boundary  are  the  forces  of  constraint  so  that  the  second 
integral  of  (2.15)  will  vanish  and  we  have 


%  -  -    jc^(TC,i^,t)UrCt,^,;t)  c/5 


(2.17) 


The  retention  of  the  second  integral  in  (2.15)  would  lead  to  the  boundary 
conditions  but  since  these  have  been  established  in  other  papers  (25,  26) 
we  will  not  derive  thera  here.   They  will  be  briefly  discussed  at  the  end 
of  the  chapter. 

Returning  to  equation  (2.12),  the  strain  energy,  with  the  help  of 
equations  (2.8),  can  be  written  as 


(r,,6,,  -^(T^y^yy  +  S^-^i 


S  -^'t- 


ax." 


^^1m  _2(r  z  ^''^- 


^i^^[{W^(§ff] 


az-dis 


(2.18) 
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Integration  over   z  yields 


H.  ^t/  /|M.te.,  +  M,,e,,  +  N,,e,, 


|N^,t,^  +  N^^t^^  +  N,^b^3 


(2.19) 


a^-ur 


^^xxr 


-M„|^.  -M^^I^.-ZM^^I^- 


lS 


where 


,  »^/i 


hll 


-h/L 


h/2. 


(2.20) 


N^3=/S^^ 


and  we  have  assumed 


K/l 


The  kinetic  energy  of  an  elastic  plate  is  given  by 


T=^|/((i%tr%^^)dl^ 


(2.21) 


■17- 


or  after  substitution  of  (2.6)  and  integration  over  z 


^4 


h|  u  t  ^  +  ^  1  + 


ds 


(2.22) 


where  /O  is  the  mass  density  of  the  plate  and  the  dot  indicates 
differentiation  with  respect  to  time.   The  second  contribution  to  the 
integral  is  due  to  rotatary  inertia.   Since  this  only  effects  the  higher 
modes  of  vibration  we  will  disregard  it.   Thus  we  write 


T 


^^ff(a'^{r'-,^')is 


(2.23) 


Combining  (2.17),  (2.19)  and  (2.23)  with  (2.10)  gives 


s 


-tz 


Ky. 


^.     s 


3X-* 


(2.24) 


Application  of  the  techniques  of  the  calculus  of  variations  to  (2.24) 
would  yield  the  equations  of  motion  in  terms  of  the  stresses  and 


displacements.   The  stress-displacement  relations  would  then  yield  the 
equations  strictly  in  terms  of  the  displacements  (see  reference  (26)). 

For  reasons  which  will  become  apparent  later  we  will,  with  the  use 
of  the  stress-displacement  relations,  obtain  the  kinetic  potential  in 
terms  of  displacement  only. 

For  an  isotropic,  elastic,  material  the  stress-strain  relations 
are 


^ax.  -   ^["^x-^K.-hC-^^ 


^^i  -"iKr^^'^-^M 


^"     E 


[^.  -Z/(^.>-+'^31^  ] 


(2.25) 


9        -    2.   _  (T 


^^^  -  ^^^   ^"^ 


Solving  the  first  two  of  these  equations  for  (ri.-j.  and  (T^j.   gives 


'3^ 


^,.^^=-    C&>.x-+^S^,)   +  At   <i;, 


_  £ 


'^r  ""  THJ 


(2.26) 


r..-   ^ 


^ 


^^"T^  ^^31"*"  ^^*-^^  "^  \^U 


r   ^.: 
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Equations  (2.26)  and  the  fourth  equation  of  (2.25)  can  be 
substituted  into  equations  (2.20)  and  with  the  aid  of  (2.8)  the  required 
integration  can  be  performed.   The  result  is  the  force-displacement 
relations 


Here, 


N    - 


_  EKC^^+^^h) 


I  +  ^ 


^4    EKCV^) 
3     ^^.  i;i 


N^   JEK_  c 

^  \'fU  3  (2.27) 


d^ur 


M,^-^{m^'W 


a'-w 


M^=  -DL%f^vf^) 


^•-r^^'-'^^ 


D=    ^^' 


\zC\-o'-) 
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We  have  a^ain  employed  the  assumption 


.h/2. 


and  further 


h/2 


-h/z. 


We  can  now  combine  (2.11),  (2.17),  (2.19),  (2.23),  and  (2.27)  to 
finally  obtain  the  kinetic  potential  strictly  in  terms  of  displacements. 
Thus, 


K- 


-^[uf  4  \r^  +  u;'    4a ur 


-eK_[  et+^J  ■>-2i^^^ejj  +  ^(\-i/)6y 


(2.28) 


ds 


Here,  the  first  term  represents  the  kinetic  energy,  the  second  the 
potential  energy  of  the  lateral  load,  the  third  the  membrane  energy  and 
the  fourth  the  bending  energy. 
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Having  the  kinetic  potential  we  could  employ  Hamilton's  principle 
and  obtain  the  three  equations  of  motion  governing  u,  v,  w.   In  the 
ensuing  analysis,  we  will  be  more  interested  in  equation  (2.28)  than 
in  the  equations  of  motion  so  that  the  latter  will  not  be  derived  here. 

The  simplification  of  equation  (2.28)  to  the  beam  equations  can 
be  made  by  setting  'VT  ,  J/*   ,  ^l^'d   equal  to  zero  so  that  we  have  for  the 
kinetic  potential  of  a  beam  undergoing  large  deformations  the  expression 


K  = 


■^.(uVur^j  n-^-^^^?f^' 


diL 


(2.29) 


This  equation  can  be  simplified  further  by  disregarding  the 
longitudinal  inertia  and  by  writing 


6      =^^ 

^  L 


L 


so    that  we   have 


K- 


PA   iJlj^  4-qixr-ii^ 


d-iipH^ 


(2.30) 


for  the  kinetic  potential  of  a  beam  undergoing  moderately  large  vibrations. 


■22- 


To  obtain  the  linearized  equations  of  motion  of  a  plate  we  disregard 
the  membrane  stresses  and  the  longitudinal  inertias  so  that  the  kinetic 
potential  becomes 


£!i^^n^-Tie^^t^)'-^^'-'^^^-%^*^^^^^] 


ds 


(2.31) 


Application  of  Hamilton's  principle  yields  (27) 

(2.32) 


3)v^ur-fpK|;ji^-.(3|^  =  c^ 


for  the  equation  governing  the  small  vibrations  of  an  elastic  plate. 
A  linear  viscous  damping  term  S^^   has  been  added  to  account  for  the 
damping  phenomena. 

2.3.   Boundary  Conditions 

For  equation  (2.32),  governing  the  small  vibrations  of  elastic 
plates,  two  conditions  are  needed  on  w  along  the  entire  boundary. 
For  example,  the  boundary  conditions  at  a  clamped  edge  are  w  =  0  and 
dlu^^n.   =  0  where  n  is  the  direction  normal  to  the  cylindrical  boundary 
in  the  xy  plane.   A  simply  supported  edge  has  the  conditions  x^j  =   Q   and 
o  w/<)n^  =   0  along  the  boundary.   Other  boundary  conditions  are  given  in 
(27). 
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When  we  consider  the  large  deflections  of  plates  as  governed  by 
equation  (2.24),  additional  boundary  conditions  must  be  prescribed.   In 
addition  to  the  two  conditions  on  w  we  need  one  condition  on  u  and  one 
condition  on  v  along  the  entire  boundary.   The  case  of  a  clamped  plate 
is  the  simplest  for  then  we  have  u  =  v  =  0.   For  a  simply  supported 
edge  we  have  Uj^  =  0,  where  u^  is  the  longitudinal  displacement  normal 
to  the  cylindrical  boundary  in  the  xy  plane.   The  other  boundary  condition 
is  obtained  from  the  condition  that  the  shear  stress  along  the  boundary 
cr^<^   vanishes.   From  (2.25),  (2.8),  and  (2.9)  we  have 


^"■bcr 


(2.33) 


Since  u   =  w  =  0  along  s  then 


(T   =  C-^^^^-   c3Ur^ 


(2.34) 


r»% 


2L         ^. 


along  s  and  the  condition  on  u   at  the  boundary  is  <3U^/<)n  =  0. 

Other  types  of  boundary  conditions  can  be  considered.   For  example 
see  (25). 


CHAPTER  III 

THE  RESPONSE  OF  LINEAR  SYSTEMS  TO  RANDOM  EXCITATION 

The  basic  aim  in  this  chapter  is  to  review  a  method  of  solution 
of  equations  of  the  type 

(3.1) 
L(ur)  -h  p  lAT  -+-  YYVW-  -Q^c->^,.*^) 

where  L  is  a  linear,  spatial,  differential  operator,  the  dot  indicates 
differentiation  with  respect  to  time  t,  and  q(x,t)  is  a  function  which 
is  not  completely  deterministic. 

We  first  must  specify  what  we  mean  by  solving  the  equations.   Since 
the  load  q  is  not  deterministic,  i.e.,  only  certain  statistical  properties 
of  it  are  known,  then  quite  naturally  all  we  can  expect  to  know  of  w  are 
certain  statistical  properties.   The  process  or  experiment  for  which  w 
is  the  result  is  said  to  be  a  stochastic  process  and  w  itself  is  called 
a  random  variable. 

Before  attempting  to  solve  equation  (3.1)  it  seems  natural  to  first 
give  a  brief  introduction  to  the  mathematical  descriptions  of  stochastic 
processes.   More  complete  and  rigorous  descriptions  of  such  processes 
can  be  found  in  any  number  of  books  (28,  29,  30). 

3.1.   Stochastic  Processes  and  Probability  Theory 

From  the  mathematical  point  of  view  a  stochastic  or  random  process 

12  N 

is  a  collection  of  functions  y  (t),  y  (t),  .  .  .  y  (t)  for  which  there 
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exists  a  probability  measure.   Each  y''"(t)  is  called  a  sample  function  or 
record  and  the  entire  collection  is  called  an  ensemble.   We  relate  this 
mathematical  model  to  the  following  physical  model.   Consider  a  certain 
experiment  which  can  be  repeated  under  similar  conditions  a  large 
number  of  times,  e.g.,  the  thermal  noise  arising  across  a  set  of  identical 
resistors.   The  outcome  of  each  experiment  is  a  different  function  y  (t) 
so  that  we  say  the  function  y(t)  which  we  are  measuring  is  a  random 
variable.    On  any  given  trial  we  cannot  predict  the  outcome  so  that  only 
certain  statistical  information  concerning  the  process  y(t)  can  be 
determined . 

The  basic  functions  that  define  a  random  process  are  the  following 
set  of  probability  density  functions: 

^l(yi''-l)  '^yi  ^  probability  of  finding  y  in  the  range  (y^^ ,  y-^   +   dy^^)  at 
t  ime  t-^. 

^2(yi>ti;  y2''^2)  '^yi  "^y2  ^  joint  probability  of  finding  y  in  the  range 
(yi>  yi  +  <^yi)  ^^    time  t-|^  and  in  the  range  (y2,  y2  +  "^y?)  ^^  time  t2. 
W3(yi,t]^;  y9,t2;  Y^'^^^    ^^l    '^yi    "^y^  ^  joint  probability  of  finding  y  in 
the  range  (yi ,  yi  4-  dyj^)  at  time  t,,  in  the  range  (y2,  y2  +  ^72^    ^^    '-i™^ 
t2  and  in  the  range  (y3,  y-^  +  dyo)  at  time  t^. 

We  continue  on  in  this  way  indefinitely.   These  set  of  functions 
must  fulfill  the  following  conditions. 

1.    Wn.^0 


/^^.•••/WnC^.,t,,...^.,:toa.^.  -i. 


The   use   of   the   word  variable    for   these   types   of    functions    is 
traditional . 


•26- 


3 .  go         oo 


W„vCl^.>!3...y,rx.^^J=Jci(^^^,  •••jWn(^,,i,;...^«..-b^;...^„.t^)dL^^ 


-oo 


since  each  W^^  implies  all  previous  \iJ^. 

It  may  happen  that  under  a  shift  of  the  t  axis  the  functions  are 
unaffected.   Such  a  process  is  called  stationary  and  we  have  Wj^Cy-j^)  dyj^  = 
probability  of  finding  y  in  the  range  (y^^ ,  yj^  +  dy^)  .   W2(y]^,y2;  ^2'^0 
dy]^  dy2  =  joint  probability  of  finding  a  pair  of  values  separated  by  a 
time  t2  -  t]^  in  the  range  (yj^,  Ji   +   dy^)  and  (y2,  y2  +  ^72)-   An"^  so  on. 
For  experimental  work  the  condition  of  stationarity  is  almost  a  necessity. 
It  is  tantamount  to  assuming  a  steady  state  condition,  i.e.,  all  transients 
of  the  system  have  disappeared. 

Quite  often  it  is  necessary  to  deal  with  more  than  one  random 
process.   That  is,  we  may  be  concerned  with  several  random  variables 
y,(t),  y2(t),  .  .  .  y^Ct).   Defining  the  process  we  then  have  the 
following  probability  density  functions: 

"l(yil.  YZI'  •  •  •  ^Nl^  tj^)  dy^^  dy2i  .  .  .  dy^^^  =  probability  that 
y^^  falls  in  the  range  (yn,  Jn   +   dyj^j^),  y2  in  the  range  (y2i>  721  "*"  ^^^21^' 
.  .  .  yj^  in  the  range  (yNl,  Y^l  +  '^yNl^  ^^  ^^^^    ^T 

W2  (yii.  721'  •  •  •  yNi'  ^i;  yw  ^22-  •  •  •  ym'  ^2^  ^^^u  •  •  •  '^^ni  "^^12 

•  •  •  '^yN2  ^  joint  probability  that  y,  falls  in  the  range  (Yh,    7^^  +  "^^11^ 
.  .  .  yj^  in  tiie  range  (yj^^,  yj^^  +  "^7^1)  ^11  at  time  tj^  and  that  y^^  falls 
in  the  range  (yj^j^,  y^j^  +  dy^^j^),  ...  y^^  in  the  range  (y^^'    y^Z  '^   '^'^m'^ 
all  at  time  t2. 

And  so  on  again.   For  simplicity  we  can  use  vectorial  notation 
and  treat  the  N  variables  as  components  of  an  N  dimensional 
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vector.   Then  in  place  of  the  above  we  may  write  W,  (y,  ,t),  W2(y]^,t,; 
y^,t2),  etc.,  with  Tj^  =  (y^^^,  y^^,  .  .  .  yNm)  • 

One  of  the  commonest  and  most  useful  density  functions  is  the 
normal  or  Gaussian  function.   For  a  random  variable  y(t)  the  Gaussian 
probability  density  function  is  defined  as  (31), 


(3.2) 


n-l    m-l 


UirrMAl 


'/i 


where 

and  where  ./V  is  the  covariance  matrix  of  elements 

The  symbol  ^   ^  indicates  ensemble  average  and  R(tj^,t  )  is  the 
correlation  function  defined  by 

If  the  process  is  stationary  then 
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and  the  process  is  completely  defined  by  the  correlation  function  B.(^ ) 
and  the  mean  m.  Furthermore,  when  ensemble  averages  may  be  replaced  by 
time  averages  (known  as  the  ergodicity  property)  then 


(3.4) 


The  Gaussian  density  function  can  be  extended  to  cover  two  or  more 
stochastic  processes.   A  discussion  of  this  situation  can  be  found  in 
(31). 

When  a  constant  parameter  linear  system  is  driven  by  a  Gaussian 
random  process  then  the  output  of  the  system  is  also  a  Gaussian  random 
process.   For  such  a  situation,  knowledge  of  the  mean  and  correlation 
function  of  the  output  permits  us  to  write  down  any  multi-variate  density 
function  for  the  output  process.   If  the  input  is  non-Gaussian  then,  in 
general,  the  output  is  non-Gaussian  and  the  mean  and  correlation  function 
no  longer  completely  define  the  process.   For  such  a  situation  no  general 
method  exists  for  finding  the  probability  density  function  of  the  output. 

One  other  useful  function  in  stochastic  processes  is  the  power 
spectral  density.   It  can  be  defined  as  the  Fourier  transform  of  the 
correlation  function.   Thus, 


,U)V  (3.5) 


^  Lu)=Jv((z]e^'^  dz 
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The  quantity  F(a))  dcO  is  the  amount  of  power  in  the  frequency  range 
((j),  oJ  +  du>)  and  hence  the  name  power  spectral  density. 
The  inverse  transform  is 


atJ 


(3.6) 


Note    that 


"R  (o)  -  i-„-/°FC<j)  Acj  =  I^Tooiri^'^'^  ^ 


(3.7) 


=  total  power  of  the  process. 

3.2.   Response  of  Linear  Systems 

Since  the  correlation  function  for  linear  Gaussian  systems  is  so 
important  we  will  outline  here  a  method  of  obtaining  this  function. 
This  will  later  permit  comparison  between  linear  and  nonlinear  theories 
for  mean  squared  values  of  displacement. 

We  start  with  the  equation  governing  the  system  in  the  form  of 
(3.1).   We  formally  seek  a  solution  in  terms  of  normal  modes.   That  is, 
we  consider 

LC^JJ")  +  7A\JJ-  -=   O  ^^-^^ 

and  seek  solutions  satisfying  the  boundary  conditions.   Let  these  be 


J*-^**^^  C3  9') 
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whereUA   is  the  natural  frequency  of  the  n  ^  mode.   Substituting  (3.9) 


'On 
into  (3.8)  yields 


L  (oCj  -  CJol  mUr^  (3.10) 


We  formally  seek  a  solution  to  (3.1)  in  the  form 


(3.11) 


YV=» 


Substitution  of  this  into  (3.1)  gives 


1  t  •  •  • 


) 


(3.12) 


Invoking  (3.10)  reduces  this  to 


(3.13) 


mX^c^.f^-Hp^o^4^  +  m^  U^c^^^^^  c^C^,-i) 


A  basic  property  of  normal  modes  is  that  of  orthogonality, 


flino(^C%.1o(^CX)cix  -  (S^^^VvA^  (3.14) 


where  dwir)  is  the  Kronecker  delta. 

Thus,  multiplying  (3.13)  by  o(^(x)  and  integrating  over  all  x 
gives  the  result 
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^^+^?^+CJo^^,^-^Jc^C^.-t)c<.C^)cU-p.Ci)    (3.15) 


The  solution  to  (3.15)  may  be  written  in  the  form 


where 


^■^^''^^/'^iS^Tzzr'^^ 


-m 


the  inverse  Fourier  transform  of 

The  function  ^^(T)    is  the  impulse  response  function  of  equation  (3.15) 

From  (3.11),  (3.15),  and  (3.16),  we  have,  as  a  formal  solution 
for  each  sample  response  w(j<,t),  the  result 


(3.18) 


where  the  lower  limit  on  h  ('2")  has  been  changed  to  -co    since  this 
function  is  zero  for  2"^  0- 
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The  cross  correlation  function  is  the  statistical  average  of 
w(x,t)  and  w(x',t  +t )  which  from  (3.18)  is 

CO 

-oo 

LL 

The  function 

h^i\.^.^.:^^^--<\Ln,:i:?\\L'^:^.^^  (3.20) 

is  the  cross  correlation  function  of  the  load  q  and  it  completely 
determines  the  cross  correlation  of  the  response  "Uj  via  equation  (3.20). 
An  important  special  case  of  CD  is 

^'K.^.'^.Al  =  ^  ^l\-^)  ^(^.-rj  (3.21) 

where  o  is  the  Dirac  delta  function.  This  permits  no  correlation  in 
space  nor  time  and  is  obviously  not  physically  realizeable.  However, 
for  lightly  damped  systems  it  represents  a  fair  approximation  to  reality. 
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Substitution  of    (3.21)    into    (3.19)    gives 


-oe 

For  the  important  case  of  lightly  damped  systems  we  have 
for  all  n  so  that  (3.17)  becomes 


^     ^   ^   XU^  gJo^^  -^-^^  (3.23) 


Substitution  of  (3.23)  into  (3.22)  would  yield  the  cross  correlation 
function  of  the  response.   Upon  setting  x  =  x' ,  "Z"  =  0,  we  obtain  the 
mean  squared  response  at  the  point  x 

9°   ,   .2.  /  .\  (3.24) 
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For  a  simply  supported  beam. 


IvVroc 


r^^  -^ 


5L 


so    that 


i    "° 


(3.25) 


For  systems  with  more  than  one  space  variable,  e.g.,  plates  or 
shells,  the  method  of  attack  is  essentially  the  same.   It  hinges  on  the 
method  of  normal  modes.   For  complicated  structures  these  functions  are 
sometimes  difficult  to  obtain  and  approximate  techniques  must  be  employed, 


CHAPTER  IV 


THE  FOKKER-PLANCK  EQUATION  AND  ITS  APPLICATION  TO  SOME 
NONLINEAR  LUMPED  PARAMETER  SYSTEMS 


The  classification  of  random  processes  is  discussed  in  this 
chapter.   It  is  shown  how  the  assumption  of  a  Markoff  process  permits 
derivation  of  the  Fokker-Planck  equation  governing  the  conditional 
probability  density  function  of  a  set  of  stochastic  variables.   It  is 
then  shown  how  this  equation  has  been  used  to  obtain  the  stationary 
first-order  probability  density  function  governing  some  nonlinear 
lumped  parameter  systems. 

4.1.   Classification  of  Random  Processes 

In  order  to  discuss  the  classification  of  random  processes  it  is 
necessary  to  introduce  the  concept  of  conditional  probability  density 
functions.   We  define  these  functions  in  the  following  manner  (see 
reference  31) : 

^2  (yi'^^l  I  T2''^2)  "^72  ~  probability  that,  if  y  has  the  value  y^^  at 
time  t,  then'y'will  have  values  in  the  interval  (y'2,  72  +  ^¥2^    ^^    '■^'"^ 

^2  (^2^y   t^l)- 

Pn  (n-'^i;  y2'^2;  •  •  •  %-l'    t^-i  1  %,    t^)  ^y^,  =  probability  that  if 
"y  has  the  values  "y^,   *y2»  •  •  •  Tn-l  ^^  the  respective  times  t^,  t2,  .  . 
^n-1  then  "y^  will  have  values  in  the  interval  (y^,    y^  +   dy^^)  at  time 
tn  (tn^  t^.^,  .  .  .-^^t^). 
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With  these  definitions  we  will  then  have 

and  so  on. 

The  Pn  must  fulfill  the  following  conditions 


(4.1) 


09    ^ 


(4.2) 

which  follow  from  definition. 

Some  further  important  properties  of  the  conditional  probability 
density  function  are: 

-t:^o''^^  ^^-'^''i^.*'') = S(^.-^,)  (4.3) 
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so  that 


(4.4) 


since  it  is  certain  that  yo  =  y]^  at  t2  =  t^^.   Here  o  (y2  ~  y^)  is  the 
Dirac  delta  function. 
We  also  have 

so  that 

since  there  is  no  statistical  dependence  between  values  of  y  observed 
at  times  sufficiently  separate. 

We  are  now  ready  to  start  classification  of  random  processes. 
The  simplest  type  of  process  is  the  purely  random  process  for  which 
we  have 


'Pn.(^>,T^.i^.,i.-...'i^>,..^-)c,.A^-v,in^-^,C^rv,^^^ 


(4.7) 


so    that    from    (4.1)    we   have 


(4.8) 


■\/\i, 
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This  last  equation  tells  us  that  the  purely  random  process  is  one  in 
which  any  y^^  and  y   for  t^^  ?^  t^^  are  statistically  independent. 

The  next  more  complicated  process  is  known  as  a  Markoff  process. 
It  is  the  situation  in  which  all  the  information  is  contained  in  the 
second-order  probability  density  function  Wj  (yi)ti^  y9,t2).   For  the 
definition  of  the  Markoff  process  we  have 

This  equation  tells  us  that  the  probability  that  y  has  values  in  the 
range  (yn,  Yn  +  ^^Yn^  ^^  time  t   given  that  it  takes  on  the  values  yj^, 
y2,  .  .  .  Yn-l  ^^  times  tj^,  t?,  .  .  .  t^^.j^  respectively  depends  only 
on  the  value  of  "y  at  the  previous  time  t^_i. 
Substitution  of  (4.9)  into  (4.1)  gives 


(4.10) 


^nC^s^.:^.A.V-^n,i.VNKX^.^MR(M^t-.WAfw.-it^ 


so  that  the  process  is  completely  specified  by  P^  (yn-1  >  '■n-l  I  ^n'  ^n^ 
since  Wj^  (YIi^i)  is  found  from  the  relation 


_  (4.11) 

1. 

We  can  continue  on  in  this  way  for  more  complex  problems.   Thus 

^3  (Yn-2>  J^n-2;  Yn-1'  '^n-l  \   "Yn'  '^n^  defines  the  next  more  complicated 
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process,  P^  (Yn-S'  "^n-B'  ^"n-2'  ^n-ZJ  Tn-l^  ^n-l  \    Tn-  ^n)    the  next,  etc. 
However,  in  this  analysis  it  is  the  Markoff  process  with  which  we  are 
concerned.   Therefore  the  higher-order  processes  will  not  be  discussed. 

It  might  seem  that  equations  (4.2)  are  the  only  restrictions  on 
^2  (yi'  ^l   \     y?'  ^2).   However,  for  a  Markoff  process  this  is  not  the 
case.   P2  (yi,  ^l   \     y2 >  ^2)  must  satisfy  the  Smoluchowski  equation, 
which  we  will  now  establish  (see  reference  31) .   We  start  from  the 
equation, 

integrate  over  y^  and  employ  (4.9)  so  that  we  have 


But, 


VJ 


i(^,  ,i,  ;"i|.,i.>'^>  C^-,''-^'^'-  C^m-^A^m"^-"^  (*•"'> 


so  that  upon  using  this  in  (4.13)  we  arrive  at  Smoluchowski' s  equation 
in  the  form 


Ti  c^,.i,i^.,-i.')  --  J^x\M'^^^o)\q.^>\^.^.)  i^. 


(4.15) 
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This  is  the  basic  equation  of  a  Markoff  process.   In  the  next 
section  we  show  how,  with  the  proper  assumptions,  this  equation  can  be 
used  to  derive  the  Fokker-Planck  diffusion  equation. 

4.2.   The  Fokker-Planck  Equation 

Before  actually  deriving  the  Fokker-Planck  equation,  let  us  look 
at  the  differential  equations  with  which  we  will  ultimately  be  dealing. 
Consider  an  N  degree  of  freedom  system  governed  by  the  differential 
equation 

where  y  =  (yi,  y2,  •  •  •  YN'  yi>  VZ'  •  •  •  Yn^  '  ^m  ^^^  ^^^^   ^   variables 
with  which  we  are  concerned  and  y^^^  are  their  derivatives  with  respect  to 
t.   Further  F  is  a  deterministic  vector  valued  function.   This  equation 
would  produce  a  deterministic  trajectory  in  the  2N-dimensional  phase 
space  and  would  be  completely  determined  by  equation  (4.16)  and  the 
specification  of  y  at  some  time  t^ . 

Let  us  add  to  equation  (4.16)  a  stochastic  forcing  function  f(t) 
so  that  we  now  have  a  set  of  equations  governing  each  sample  of  y  which 
have  the  form 

obt  ^ 

The  trajectory  of  the  phase  space  is  now  a  stochastic  process.   Clearly, 
the  position  y(t2)  at  the  end  of  any  infinitesimal  interval  of  time 
depends  only  on  the  value 'y'(t^)  at  the  beginning  of  the  interval  (t2,t,) 


-Al- 


and on  the  stochastic  forcing  function  f  acting  during  this  interval. 
We  now  make  the  assumptions  that  the  forcing  function  f  is  Gaussian 
with  zero  mean.   We  further  assume  that  the  forcing  functions  acting 
on  the  system  at  any  two  small  consecutive  time  intervals  are 
statistically  independent.   These  assumptions  make  it  necessary  for  us 
to  take 

<L.  Ci^>  =0  (4.18) 


<\rr.ih)   |.  a+^^  >  -  l?,n.  ^Ct)  (4.19) 

where  fj^  is  the  ra"^^  component  of  f  and  I^^  ^^  some  function  of  m  and  n. 
Furthermore,  since  the  position  at  the  end  of  an  infinitesimal  time 
interval  depends  only  on  the  value  y(tj^)  at  the  beginning  of  the  interval, 
and  on  the  forcing  function  acting  during  the  interval,  which  according 
to  our  assumption  is  independent  of  the  forcing  function  acting  outside 
the  interval,  then  the  trajectory  of  the  phase  space  is  a  Markoff 
process.   It  is  completely  defined  by  the  conditional  probability  density 
function  P2  (yi.t]^  |  y7,t2),  which  must  satisfy  the  Smoluchowski 
equation  (4.15). 

With  the  proper  assumptions  equation  (4.15)  can  be  used  to  derive 
the  Fokker-Planck  equation  governing  P2  (^l,t]^  \  y2,t2).   The  derivation 
can  be  found  in  many  places  (4,  10,  31).   Here  we  follow  the  derivation 
given  in  (10) . 

To  start  the  derivation  of  the  Fokker-Planck  equation  we  consider 
the  first  and  second  moments  of  the  displacement  of  the  phase  point  in 
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an  infinitesimal  time.   These  are 

We  assume  that  these  are  of  order  A  t  and  that  all  higher  moments  are 
of  higher  order  of  ^t.   The  first  assumption  insures  the  existence  of 
the  following  limits. 


(4.21) 


It  has  been  pointed  out  in  (4,  10)  that  these  assumptions  are  tantamount 
to  the  assumption  of  a  Gaussian  process  for  the  disturbances. 

Having  made  these  preliminary  assumptions  we  consider  an  arbitrary 
scalar  function  R(y)  which  vanishes  sufficiently  fast  to  zero  at  infinity. 
Multiplying  the  Smoluchowski  equation  (4.15)  by  this  function  and 
integrating  over  the  entire  phase  space  gives 

(4.22) 

and  we  have  interchanged  the  order  of  integration.   We  now  develop  R(y) 

-»  .^ 
in  a  Taylor  series  in  (y  -  x) , 
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Upon  using  (4.23),  (4.21)  and  the  assumptions  concerning  the  higher 
moments,  the  right-hand  side  of  (4.22)  becomes 


JRc^)T.c^i\^.,t«)<u 


(4.24) 


VAi: J  [||  \  i. .  CU )  |^3?A^.  .^o\^ ,i^  c)J^ 


-¥    o 


A-k: 


Integrating  by  parts,  writing  y  for  x  and  putting  the  result  in 
(4.22)  gives 


(4.25) 


ZW  Zfi 


+Tll._,-^::^„l^-^?Ai,.t«\a,i)lld| 


\  OA^ 
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Taking  limits  as  A  t— •  0  yields 


i^^ii^^tk''''^'^'^ii.^^^^-^%--^  <--> 


Since  R(y)  is  arbitrary,  the  bracketed  expression  must  vanish, 
leaving  us  with  the  Fokker-Planck  equation. 

This  is  a  parabolic  diffusion  equation.   The  required  solution  is 
the  positive  one  with 

If  all  transients  of  the  system  have  died  out  and  a  steady  state 
condition  has  been  reached,  then  Piyo^t^^  y,t)  — *  Wj^(y,t)  so  that  the 
Fokker-Planck  equation  becomes 

It  still  remains  to  be  shown  how  equation  (4.29)  can  be  used  to 
determine  W^^  governing  the  variables  of  equation  (4.17).   The  connection, 
of  course,  is  through  the  moments  d   and  d   .   Indeed  we  had, 

JL^C^A^-  ^^-^^$C^^-M^jivlC^,\^X,i-^rd.^cU  (4.30) 
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Now  at  the  beginning  of  the  time  interval  x^^  had  the  value  y^  and  at 
the  end  the  value  y^,  +  ^  Ym  ^°  that  dx^  =  d(/:iyj^)  and  we  have 

^  Ic/^  _L_-r  ^^^^^^/^(^£,^  j/^l'jcJLCAp   (4.31) 

Similarly,  it  can  be  demonstrated  that 

dLmT^   -  U^     _L_<AM^^M^>  (^-32) 

These  moments  can  now  be  computed  from  the  differential  equation 
(4.17)  and  thus  the  corresponding  Fokker-Planck  equation  is  fully 
derived.   The  method  is  best  illustrated  by  examples.   Therefore,  to 
close  this  chapter  we  consider  two  lumped  parameter  systems  which 
have  been  analyzed  by  this  method. 

The  first  system  is  a  one-degree  of  freedom,  nonlinear  oscillator. 
The  equations  of  motion  are 


^•+^.t,.^).|(t) 


(4.33) 


Also 
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Writing  y^  =  y,  y2  ~  y  equation  (4.33)  is  equivalent  to  the  two  equations 


The  coefficients  of  the  Fokker-Planck  equation  are  therefore 


A-t^o 


(4.35) 


'^-^V^,ii<n^^^-^^--''^1>^ 


(;t  11  -  Xwn. 


=  0 


(4.36) 


Substitution  of  these  coefficients  into  equation  (4.29)  yields 


The  solution  to  this  equation  as  given  in  (9)  is 


(4.38) 
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where  C  is  a  normalizing  constant.   If  k(y)  is  linear,  then  this 
reduces  to  the  Gaussian  density  function  as  it  should. 

The  second  illustration  of  the  application  of  the  Fokker-Planck 
equation  to  lumped  parameter  systems  consists  of  a  loaded  nonlinear 
string  as  analyzed  by  Ariaratnam  (11).   The  N  equations  governing  the 
system  are  of  the  form 


(4.39) 

To  simplify  these  equations,  the  deriving  functions  and  the  response 

are  expanded  in  terms  of  the  eigenfunctions  of  the  linear  problem.   Thus, 


^-^^ 


(4.40) 


\ 


Substitution  of  (4.40)  into  (4.39)  yields 


m. 


(4.41) 


This  equation  may  be  used  to  derive  the  moments  appearing  in  the  Fokker- 
Planck  equation.   The  solution  of  the  resulting  equation,  as  obtained 
by  Ariaratnam  is 
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.  N  \^ 


(4.42) 


where  C  is  a  normalizing  constant  and 

This  expression  can  be  used  to  obtain  the  mean  squared  displacement  of 
the  various  masses  on  the  strong. 

The  Fokker-Planck  equation  has  been  used  in  the  past  to  solve 
several  nonlinear  lumped  parameter  systems.   In  the  next  chapter  the 
responses  of  some  nonlinear  continuous  structures  such  as  beams  and 
plates  are  investigated  by  this  method. 


CHAPTER  V 


APPLICATION  OF  THE  FOKKER-PLANCK  EQUATION  TO  NONLINEAR 

ELASTIC  SYSTEMS 


In  this  chapter  it  is  shown  how  the  Fokker-Planck  equation  can  be 
used  to  investigate  the  finite  responses  of  plates  which  have  been 
subjected  to  white  noise  excitation.   A  general  solution  to  the  Fokker- 
Planck  equation  is  given  which  is  applicable  to  plates  with  any 
boundary  conditions.   Detailed  solutions  are  presented  for  a  simply 
supported  beam  and  a  simply  supported  plate. 

5.1.   General  Theory 

In  Chapter  II  it  was  shown  that  the  forced  vibrations  of  an  elastic 
plate  are  governed  by  the  following  equation 

S  j  Kdi     -  0  (5.1) 

where 


K=  f[[^(Al^+vNw^U^w-\/,ys 


(5.2) 


and  Vg  is  the  strain  energy. 

Instead  of  applying  variational  techniques  to  equations  (5.1)  and 
(5.2)  and  thus  obtaining  the  three  equations  of  motion  governing  u,  v,  w 
we  proceed  as  follows:   We  expand  each  sample  function  of  w  and  q  in  a 

-49- 
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series  of  the  eigenfunctions  of  the  linear  problem.   This  is  valid  as 
long  as  each  sample  of  w  and  q  has  continuous  derivitives  up  to  fourth 
order  (reference  32,  p.  370),  a  condition  we  now  assume.   Thus,  we  write 


N  M 


%^'-'^'^^-l,lA-^^'''^-^'--'i^ 


(5.3) 


N  Kl 


<1     -n=i  m-l  *J 


(5.4) 


The  infinite  series  has  been  terminated  at  some  N,  which  is  later  to  be 
specified.   Of  course  this  invalidates  the  equality  sign  of  equations 
(5.3)  and  (5.4).   However,  as  long  as  the  infinite  series  representing 
q  and  w  converge,  the  finite  sum  can  be  made  as  accurate  as  desired  by 
properly  selecting  N. 

Similar  to  equation  (3.10)  the  eigenfunctions  of  the  plate  must 
satisfy  the  equation 


L^*'  X^„.W^.  =  0 


(5.5) 


where  the   X   are  the  eigenvalues  determined  from  the  frequency 

'  mn  '^  ^ 

equation.   In  addition  the  W   must  satisfy  the  appropriate  boundary 
conditions. 

Each  sample  function  of  u  and  v  is  also  expanded  in  an  infinite 
set  of  functions.   Wc  choose  this  set  to  be  orthogonal,  to  satisfy  the 
boundary  conditions  and  to  be  such  that  term  by  term  differentiation 
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of  the  infinite  series  is  possible.   Thus, 


(5.6) 


These  infinite  series  have  also  been  terminated  with  the  previous  argument 
concerning  convergence  still  applying.   For  convenience  the  same  N  has 
been  chosen  for  all  series. 

An  example  of  a  proper  expansion  would  be 


\JlU-»VJ,-L^   "=^11   JJLtAnCi^    AW»  SilLt     /iAM 


tOT^         (5.7) 


for  a  clamped  rectangular  plate  of  sides  a  and  b  and  with  the  origin  of 

the  coordinate  system  at  a  corner  of  the  plate.   This  would  be  a  double 

,      .   m'^  X      n  <t  Y       m  IT  X    .   nTT  y 
Fourier  series  with  terms  such  as  sin  -         cos  — r-^  ,  cos  — - —  sin  — g-^, 

cos  SLllii  cos  "5  y   omitted.   These  terms  can  be  omitted  if  we  consider 
a         b 

the  extension  of  u  onto  the  intervals  (-a,0)  and  (-b,0)  to  be  an  odd 
function.   Now  since  (5.7)  is  a  Fourier  sine  series,  which  is  zero  at 
X  =  0,  X  =  a,  y  =  0,  y  =  b,  then  it  is  continuous  throughout  the  entire 
xy  plane  and  may  be  differentiated  term  by  term  (see  Theorem  I  on  page 
137  of  reference  32).   Therefore,  for  a  clamped  rectangular  plate 
equation  (5.7)  is  a  series  expansion  of  u  which  satisfies  the  boundary 
conditions  and  the  conditions  for  at  least  one  differentiation  provided 
each  sample  u  is  integrable  in  the  Lebesgue  sense,  a  condition  we  now 
assume . 
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If  we  now  substitute  equations  (5.3),  (5.4)  and  (5.6)  into  equation 
(5.2)  and  perform  the  indicated  integration,  we  obtain 


(5.8) 


-  ^y^y,   ^^-n~'^(  U..,..  ^^»,.M.j.  '^Hv*^ 


where  V  is  the  strain  energy  V^  after  integration  and 


^'-'^  ~  j[V^\(?CU)iS  (5.9) 


Here  the  orthogonality  property  of  the  eigenfunctions  has  been  used. 

To  satisfy  equation  (5.1)  we  consider  the  u   ,  v   ,  w   as 
generalized  coordinates  and  apply  the  Euler-Lagrange  variational 
equations,  which  are  of  the  form 


(5.10) 


where  the  9^   are  the  generalized  coordinates.   Application  of  (5.10) 
to  (5.8)  yields  the  following  set  of  differential  equations 
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\r r— ; —  ^ —  (5.11) 


where  we  have  introduced  the  same  linear  viscous  damping  term  of  the 
linear  problem  (see  equations  (2.41)  and  (3.15)). 

We  now  write  (5.11)  as  a  set  of  first-order  equations  by  setting 

Equations  (5.11)  then  become, 


(5.12) 


XL     =  --1—    SX 


'vr     -    AT 


(5.13) 
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^    =  -  ^ ^v 

rr\T\ 


ptT-bvVsVV  "^ATy 


AJj'w.^  -   iat^ 


-im  n 


'i^-m.N- C7 —  r:T-  -  L^ — :  4--r-  9.^'^ 


These  equations  constitute  a  set  of  stochastic  differential  equations. 
They  are  of  the  same  form  as  (4.17).   The  vector  y  has  the  components 

"mn.  u^n'  %n'  ^mn'  ^mn'  ^mn  ^°^   ^^^   values  of  m  and  n.   This  is  a 

2 
total  of  6N   components.   The  stochastic  driving  functions  f   are 

equivalent  to  the  qj^j^.   Therefore,  if  we  assume  that  the  q   belong 

to  a  Gaussian  random  process  and  that 

(5.14) 

then  the  arguments  of  Chapter  IV  hold  and  the  trajectory  in  the  phase 
space  of  equations  (5.13)  constitute  a  Markoff  process.   As  was 
demonstrated  in  Chapter  IV,  this  implies  that  the  stationary  probability 
density  function  W]^(y)  must  satisfy  the  following  Fokker-Planck  equation 
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where, 


dl^-   luw^-^ 


Equations  (5.14)  imply  that 


(5.16) 


<QL^,^^-V)o^L<,^^h^'^)>~-\\Li^.^,y:,^)%>(X^         (5.n) 


7ith 


(5. IK) 


The  moments  d^  and  d    of  equation  (5.16)  can  be  calculated  from 
the  differential  equations  (5.13).   For  bookkeeping  purposes  the 
following  notation  is  introduced: 


Adt->o  At 


<  A  U,^^> 


cL__  -  Mjyy\ L_ 


■»A'Y>^ 


A-L-*C)  aI- 


<  AU_^^> 


ot 


(5.19) 
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A-L-*C>   ^*- 


n,    u. 


C^mrur-s^    U^   _1^<  ALL^,  Aa,,> 


dmtvrs     —     i^'^ 


Ai:*o   Ai. 


<.  AU,^^  /^Uvs> 


Then,  with  the  aid  of  equations  (5.13)  these  moments  can  be  calculated 
to  be 


oL  m  TV  -    Vv 


vnrv 


A.^      -    =^ ^^  (5  20) 
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ioii^I  /  v^(-n--^"^d^^> 


*^     t.  ^*-  K^ 


"^^         »^  r»\  f  rs 


all   other   d^^rs   =  0- 

Substitution  of  these  expressions  into  equation  (5.15)  yields  the 
following  Fokker-Planck  equation: 
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This  is  indeed  a  complex  equation  for  which  there  exists  no 
standard  technique  for  obtaining  the  solution.   By  a  process  of  trial 
and  error  the  author  was  fortunate  in  obtaining  a  solution  for  the  case 
in  which 


tNS 


(5.22) 


This  imposes  the  restriction  that  the  load  is  completely  uncorrelated 
in  time  and  space,  i.e., 

The  solution  to  (5.21)  with  the  restriction  of  (5.22)  is 


+b^^Ar^^  -tc, 


.TV^^^O^]! 


(5.24) 


where  C'  is  a  constant  to  be  obtained  from  the  normality  condition. 
The  general  uniqueness  theorem  to  the  Fokker-Planck  equation  given  in 
(33)  tells  us  that  this  is  the  only  solution. 

It  is  to  be  noted  that  the  distribution  of  the  velocity  variables 
is  Gaussian.   If  we  integrate  over  these  variables  from  -  CX3    to  +oO      we 
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obtain   the   probability   density    function  of   the   u      ,    v        and  w  Thus, 


ran'   mn      mn • 


V^.(u.,..a^^.,...^xT.,...uJ^.^  -  C'siyppl-^]  (5.25) 

The  constant  C  is  obtained  from  the  normality  condition. 

We  are  ultimately  interested  in  the  probability  density  function 

of  Wj^jj^  since  from  it  we  can  determine  quantities  such  as  mean  squared 

displacement,  mean  squared  stress,  etc.   To  obtain  this  we  integrate 

equation  (5.25)  over  all  u  „  and  v so  that 

^       ^     '  mn      mn 

This  is  as  far  as  can  be  gone  in  such  generality.   To  proceed  further 
we  need  to  compute  V.   Several  important  cases  are  presented  in  the 
next  section.   It  should  be  noted  in  passing  that  when  the  load  has  some 
spatial  correlations  it  would  be  advisable  to  employ  some  approximate 
technique  to  solve  equation  (5.21)  governing  the  probability  density 
function  of  the  modal  amplitudes. 

5.2.   Some  Special  Cases 

a.   Simply  Supported  Beam 

The  kinetic  potential  of  a  beam  with  moderately  large  vibrations 
is  given  by  equation  (2.29)  so  that  the  strain  energy  is 


^^     2. 


[i^M^.  ^iteH 


(5.27) 
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The  eigenf unctions  of  the  linear  problem  are  sin  2L_JS  g^,  that  we 


write 


M 


VJ^TL^i.^   -  Y-lij;A>L^  n-TT  ?c 


(5.28) 


n=  I 


The  first  four  equations  of  (5.20)  are  satisfied  identically  since 

u  =  V  =  0  and  since  the  strain  energy  has  been  expressed  independently 

o  f  u  and  v . 

The  probability  density  function  of  the  w  is  then  given  by 


W.Cu),  ,...u)^^  -  CMfi-tl,] 


(5.29) 


with  V  given  by  (5.27)  after  substitution  of  (5.28).   Thus 


W.  Cvr„.  .U3^^  .  C' (^  -^^  Y^'< -^i^^[l"^^<l 


^^^L 


(5.30) 


where 


pZT.   1>' 


^'   = 


^ 


(5.31) 


N),  ^  K 
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Here  (TL  is  the  mean  square  deflection  of  the  first  mode  of  the  linear 
problem  obtained  by  setting  V\  =  OO.  We  have  replaced  N'/L  by  N  ,  which 
is  interpreted  as  the  power  spectral  density  of  the  average  load  acting 
on  the  beam,  i.e., 


C£> 


L\'» 


■u)l: 


L. 


(5.32) 


where  q  is  now  to  be  interpreted  as  force  per  unit  length  and  we  have 
used  the  fact  that 


<<^i,^,-k^c^i.iL,X^'X:)>-^H'  ^L%-7C)IC^^ 


(5.33) 


The  static  deterministic  counterpart  of  N^  for  a  uniformly  loaded  beam 

2 
is,  of  course,  q  ,  i.e.,  the  square  of  the  constant  load  per  unit  length, 

It  is  seen  from  equation  (5.30)  that  the  nonlinearity  of  the  beam 
causes  the  probability  density  functions  of  the  modal  amplitudes  to 
become  non-Gaussian.   Furthermore  these  variables  are  no  longer 
statistically  independent. 

The  mean  squared  response  of  the  beam  is  given  by 


M  «^ 


(5.3A) 
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where 


<'lAJ•^^J■^>  =  J•••j^^^^^J'>^W.(vJ.^...uJ^^dy,..  du^,  (5.35) 


—  ao        -oo 


Because  equation  (5.30)  is  an  even  function  of  the  w^  we  will  have, 

<\J-^\A^^>  ^  C?   ,         -WN^^  (5.36) 

so  that  the  modal  amplitudes  are  linearly  independent. 

The  mean  squared  stress  in  the  beam  depends  upon  '^V.'a'^J-  •  /^  which 
is  given  by 


,  (5.37) 


9  / 

In  the  linear  problem,  <  w^^  >   is  of  the  order  1/m   so  that  the 
infinite  series  does  not  converge.   This  was  discovered  by  Eringen  (12) 
and  attributed  to  two  factors:   (a)  The  o   functions  appearing  in 
equation  (5.23)  and  (b)  The  inadequacy  of  Bernoulli-Euler  beam  theory. 
To  obtain  expressions  for  mean-squared  stresses,  Samuels  and  Eringen  (34) 
investigated  a  Timoshenko  beam  and  Crandall  and  Yildiz  (3  5)  later 
investigated  many  different  beam  models  with  different  damping  mechanisms, 
It  was  found  that  these  more  refined  theories  produced  finite  mean- 
squared  stresses. 
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Since  we  have  only  an  integral  representation  for  ^  w   ^   we 
cannot  rigorously  investigate  equation  (5.37).   However,  there  is 
really  no  reason  to  believe  that  the  introduction  of  the  membrane 
stresses  would  cause  (5.37)  to  converge.   To  investigate  the  mean- 
squared  stresses  in  the  nonlinear  problem,  it  would  be  most  desirable 

to  consider  a  more  refined  beam  model. 

2 
An  approximate  expression  for  <  w  'y     can  be  obtained  by 

substituting  the  expression 

into  (5.35)  and  performing  the  required  integration.   The  result,  after 
some  juggling,  valid  only  for  "moderately  large  deflections"  is 


<'UJ-^>   ^ 


5^ 


\- 


L[c^Vt+^i^-+^.I^^+a]- 


(5.35) 


lu^ 


b .   Simply  Supported  Plate 

We  consider  a  rectangular  plate  simply  supported  on  all  four  sides. 
For  the  linear  problem  the  eigenfunctions  are 


W^    C>-j^^  -    A>un  r^}JlA  Ax/a  uTtM 
^  CL  b 


(5.40) 
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where  a  and  b  arc  the  lengths  of  the  sides  of  the  plate  in  the  x  and  y 
directions  respectively. 

Earlier  in  the  chapter  it  was  shown  how,  for  a  clamped  plate,  a 
double  Fourier  series  could  be  constructed  for  the  functions  u  and  v, 
which  satisfied  the  boundary  conditions  and  the  conditions  for 
differentiability.   For  a  simply  supported  plate  we  define  the  extension 
of  u  on  the  interval  x  =  (-a,0)  to  be  odd  and  on  the  interval  y  =  (-b,0) 
to  be  even.   Then  we  can  write 

M.itM,iL^-=^Y-     XL^^C^^AAAvmlvxc^  nir^     (5.41) 

^  ^ZO   ^-<=>  0-  fc> 

Now  by  Theorem  I  of  reference  (32)  it  is  evident  that  this  series  can 
be  differentiated  with  respect  to  x  or  y  at  least  once.   The  summing  of 
the  series  with  n  =  0  is  necessary  since  a  cosine  expansion  must  include 
a  constant  term.   For  convenience  we  have  also  started  the  summation 
at  m  =  0. 

For  the  function  v  we  can  proceed  along  the  same  lines.   We  may 
therefore  write  for  the  three  displacements 

AA-CS  M  ,Tb^  -  Z  1"  Vm.  C-fcT)  Ax;rv.  mXYj^  COO  y^^ 

d       r\-o  m-u  OL  ta 


'^Ct^,H,^>*I.'C  w^^Ct^  i^AAirO^STA   aa^^mx_!4. 
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For  uniformity  we  have  started  all  suras  at  m  =  n  =  0  and  terminated 
them  at  m  =  n  =  N-1  so  that  the  total  number  of  terras  of  u^^^^  is  still 

The  kinetic  potential  of  a  plate  undergoing  large  vibrations  is 
given  by  equation  (2.28).   The  strain  energy  is  therefore  given  by 


(5.43) 


Vfei-^i^^V^-^^s^^-^'-''^^'^^'^^ 


where  the  £  ,  £.  ,  £   are  given  by  (2.9).   If  ^^?e  now  substitute 
(5.42)  into  (5.43),  perform  the  required  integration  and  substitute 
the  result  into  (5.25)  we  obtain,  after  some  rearranging 


(5.44) 
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where 


D^.p^.*..  ■=  T^^-^^  ^.p.Je>.^s.  +  T^sxKp^A.^, 
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and 

a 


o 

a 

F'  „r.  -f  ocyd'^:^^^^  A>t^-p:5U<^    zCwruoDi.   dx- 
'wp'^       I  Q.  a,  a 

A„pr^/c<^!^  c^£ii  cx>^r^  to^«=2i  d^ 
-'o 

b  (5.46) 
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To  obtain  the  probability  density  function  of  the  w   alone  we 

^  ■'  ■'  mn 

need  to  integrate  out  the  u   and  the  v     If  we  make  use  of  the 
relation 


CO 


_i£^L-Ca^M^\.^^]A^.   {^^(-^^ 


(5.47) 


which  is  given  on  page  64  of  reference  (36),  then 


(5.48) 


becomes,  after  some  manipulation 
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(5.49) 
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If  we  integrate  out  the  v   ,  again  making  use  of  (5.47),  we  finally 
obtain 


(5.50) 
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i'-ou  ^Ci-v)V3«.'-  _  (LjJiLiilA 


where  as  before  N^  is  the  power  spectral  density  of  the  average  load 
acting  on  the  plate,  i.e., 

J  i,  s 
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N^^  =  N'/ab 
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Again,  the  effect  of  the  nonlinearity  is  to  cause  the  modal 
amplitudes  to  become  statistically  dependent  and  non-Gaussian.   The 
mean-squared  displacement  of  the  plate  is  given  by 


cvJ^  C^,^,i-^  -  (Tp"  =  Z II 1  <  ur^.^r.>  A^  "^ 


*■  (5.51) 


where 


CO  CO 


Having  reduced  the  mean-squared  displacement  to  quadratures,  we  must 
again  stop  and  seek  some  approximation  to  finish  the  job.   This  is  left 
to  the  next  chapter. 

As  was  the  case  in  the  beam,  the  mean-squared  stresses  for  the 
linear  plate  do  not  converge  and  so  will  not  be  investigated  here. 

The  investigation  of  the  large  random  vibrations  of  a  simply 
supported  beam  and  simply  supported  plate  have  been  reduced  to  quadratures 
in  this  chapter.   We  leave  to  the  next  chapter  the  numerical  investigation 
of  specific  cases. 


CHAPTER  VI 


NUMERICAL  INVESTIGATION 


In  the  previous  chapter  the  probability  density  function  of  the 
modal  amplitudes  of  a  plate  undergoing  large  deflections  was  derived. 
Detailed  calculations  were  presented  for  a  simply  supported  beam  and 
a  simply  supported  plate.   The  mean-squared  deflection  of  these 
structures  was  reduced  to  quadratures.   For  the  beam  an  approximate 
formula  was  developed. 

In  this  chapter  we  investigate  more  closely  some  of  the  results 
of  the  previous  chapter.   For  the  beam,  the  linear,  approximate 
nonlinear  and  integral  representation  of  the  mean-squared  displacement 
are  compared  numerically  for  a  range  of  the  parameters.   For  the  plate 
the  mean-squared  displacement  of  the  first  mode,  which  is  the 
predominant  term,  is  numerically  investigated  for  different  aspect 
ratios  b/a. 

6.1.   Simply  Supported  Beam 

In  Figure  3,  the  mean-squared  deflection  at  the  center  of  the 
beam  (T^  ,  as  determined  by  the  linear  theory,  the  approximate  formula 
(5.39),  and  numerical  integration  of  equation  (5.35)  is  plotted  against 
(TcT   with  R  =  1/2.   For  the  range  of  tf^   considered  it  was  found,  as 
in  the  linear  theory,  that  sufficient  accuracy  is  obtained  with  N  =  1. 
It  was  also  found  that  the  approximate  formula  is  valid  over  a  small 
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2. 

range  of  (J^   .   As  long  as  the  difference  between  linear  and  nonlinear 

theories  is  not  greater  than  10  per  cent,  this  formula  gives  an  excellent 
estimate  of  the  true  mean-squared  deflection.   This  is  to  be  expected 
in  view  of  the  approximation  made  in  deriving  formula  (5.39). 

Also  plotted  in  Figure  3  is  the  mean-squared  deflection  determined 
by  numerical  integration  with  N  =  1  and  R  =  1/4  and  1/8.   These  curves 
indicate  that  for  lower  values  of  the  radius  of  gyration  the  curve  of 
the  nonlinear  theory  begins  to  deviate  sooner  from  the  straight  line  of 
the  linear  theory.   Tliis  is  to  be  expected  since  with  diminishing  R  the 
role  of  bending  dimishes.   Furthermore,  if  the  beam  were  rectangular 
then  h   =  12R   so  that 


S-^     =1^1 


-^  (6.1) 


Now  for  R  =  1/2,  the  nonlinear  theory  begins  to  appreciably  deviate 
from  the  linear  theory  at  f^     =0.3  or  at  ^=0.1.   For  R  =  1/4  we 
have  appreciable  deviation  at  (Tg   =  0.009  or  S=  0.11,  and  for 
R  =  1/8  at  (Tg^  =  0.004  or  <^  =  0.14.   Since,  as  is  well  known,  the 
linear  theory  of  beams  is  valid  only  for  w/h<.<  1,  these  results  are 
reasonable. 

In  Figure  4  the  mean-squared  deflection  as  determined  from  the 

linear  theory  and  from  numerical  integration  of  equation  (5.35)  is 

2 
plotted  against  larger  values  of  Q^    for  R  =  1/2.   Of  the  three  curves, 

the  uppermost  represents  numerical  integration  with  N  =  1.   The  lowest 

2 
curve  represents  <i  wi   >   as  determined  by  numerical  integration  with 

N  =  3.   The  middle  curve  represents  the  total  deflection,  i.e., 
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2  2 

^  w-|^  ^   "^  ^  "3  ^   with  N  =  3.   Two  observations  are  evident  from 

these  curves.   Firstly,  consideration  of  the  beam  as  a  one-degree  of 

freedom  system  (i.e.,  N  =  1)  is  no  longer  valid  for  such  large  values 

of  Og   .   Secondly,  as  in  the  linear  theory,  <  wj^  '^   gives  a 

fairly  close  estimate  (a  few  per  cent)  of  the  total  mean-squared 

2 
deflection.   However,  in  computing  ■^  wi   ^  it  is  now  necessary  to 

consider  the  effects  of  the  second  and  third  modes.   That  is  to  say, 

the  nonlinear  coupling  is  so  strong,  for  the  range  of  parameters 

considered  in  Figure  4,  that  the  second  and  third  modal  coefficients 

have  a  significant  effect  on  the  mean-squared  value  of  the  first  modal 

coefficient . 

It  should  be  pointed  out  that  for  a  rectangular  beam  with  R  =  1/2 

we  have  h   =  3,  so  that  with   <r~g   =  1.0,  the  highest  value  in  the 


graph  of  Figure  3,  we  have 


1 


^^ 


6_   ^  O.b  ^^-2) 


This  value  is  not  in  excess  of  the  applicability  of  the  nonlinear  theory 
considered  nor  of  the  values  of  practical  interest. 

6.2.   Simply  Supported  Plate 

The  first  mode  is  the  predominant  one  in  the  calculation  of  the 
mean-squared  displacement  of  the  plate.   Therefore,  if  in  equation  (5.50) 
we  consider  only  the  first  mode  we  have 
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and  all  other  Lmnllll  "^  ^mnllll  ^  ^  so  that  (6.3)  reduces  tc 
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where 


'  I=b/a 

This    is   of   the    form 


where   u,  and  c<    are   easily   calculable   once    xJ ,    a   and   b   are   k 
The  mean-squared   deflection  of   the   plate    is 


(6.6) 
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Unfortunately  this  integral  is  not  tabulated  and  we  must  resort  to 
numerical  methods. 

In  Figure  5  the  mean-squared  deflection  at  the  center  of  the  plate 

2 
(^       as  determined  from  numerical  integration  of  (6.7)  is  plotted 

2 
against  (C~        for  l)   =  1/3,  h  =  1,  and  two  difference  aspect  ratios 
Po 

of  X  =  1  and   j/  =  2 .   For  larger  aspect  ratios  the  linear  and  nonlinear 
theories  give  higher  mean-squared  deflections.   This  is  what  is  to  be 
expected  since  increasing  the  aspect  ratio  is  equivalent  to  moving 
apart  one  pair  of  supports  of  the  plate. 

Also  evident  from  Figure  5  is  that  the  nonlinearity  of  the  plate 
causes  a  reduction  in  the  mean-squared  deflection.   Furthermore,  for 
the  smaller  aspect  ratio  the  percentage  deviation  of  the  nonlinear 
theory  is  slightly  greater  than  for  the  higher  aspect  ratio. 

In  this  chapter  the  mean-squared  deflection  of  the  beam  and  plate 
has  been  investigated  by  numerical  integration.   The  results  indicate 
reduction  of  these  quantities  and  a  significant  coupling  of  the  modes 
for  sufficiently  large  deflections. 
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CHAPTER  VII 


CONCLUSIONS 


Tills  paper  has  been  devoted  to  the  study  of  the  random  vibrations 
of  some  nonlinear  elastic  systems.   The  equations  of  motion  of  thin 
elastic  bodies  with  large  deflections  have  been  presented  in  variational 
form.   General  truncated  series  expansions  of  the  middle  surface 
displacements  have  been  performed  and  the  Euler-Lagrange  variational 
equations  have  been  used  to  obtain  the  nonlinearly  coupled  differential 
equations  governing  the  series'  coefficients. 

The  assumption  of  uncorrelated,  Gaussian  loading  permitted  the 
phase-space  of  the  series'  coefficients  to  be  identified  with  a  Markoff 
process.   This,  in  turn,  permitted  the  derivation  of  the  Fokker-Planck 
equation  governing  the  probability  density  function  of  the  series' 
coefficients.   Tlie  moments  for  this  equation  v;cre  obtained  from  the 
governing  set  of  differential  equations. 

A  solution  to  the  Fokker-Planck  equation  for  the  case  of  spatially 
uncorrelated  white  noise  has  been  obtained.   Several  special  cases  have 
been  worked  out  in  detail.   In  particular,  the  probability  density 
functions  for  a  simply  supported  beam  and  a  simply  supported  plate  have 
been  computed.   The  mean-squared  displacements  of  these  systems  were 
reduced  to  quadratures.   For  the  case  of  the  beam,  an  approximate 
formula  was  developed. 
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Numerical  integration  of  the  exact  expression  for  the  mean-squared 
displacement  of  the  beam  has  shown  a  reduction  of  this  quantity  as 
compared  to  the  linear  theory.   Calculations  have  also  shown  a  limited 
range  of  applicability  of  the  approximate  formula.   Furthermore,  it  was 
shown  that  the  first  mode  still  gives  a  good  estimate  of  the  mean-squared 
response  but  the  coupling  effect  of  the  modes  is  so  important  for 
sufficiently  large  deflections,  that  the  effect  of  the  higher  modes  must 
be  taken  into  account  when  computing  the  mean-squared  value  of  the  first 
mode.   Numerical  integration  also  showed  a  reduction  of  the  mean-squared 
deflection  of  the  plate. 

While  this  paper  has  presented  a  method  of  attacking  the  random 
vibrations  of  elastic  systems  with  large  deflections,  it  must  be 
pointed  out  that  tlie  solutions  presented  are  valid  only  for  completely 
uncorrelated  loadings.   This  is,  of  course,  physically  unrealizeable . 
Nevertheless,  the  results  of  this  analysis  should  give  some  insight  into 
the  problems  of  the  random  vibrations  of  nonlinear  elastic  systems. 


LIST  OF  REFERENCES 


1.  Crandall,  S.  H.,  ed .   Random  Vibration.   Vol.  1,  the  M.I.T.  Press, 

Cambridge,  Massachusetts,  1948. 

2.  Crandall,  S.  H.,  ed.   Random  Vibration.   Vol.  2,  the  M.I.T,  Press, 

Cambridge,  Massachusetts,  1963. 

3.  Unlenbeck,  G.  E.,  and  Ornstein,  L.  S.   "On  the  Theory  of  Brownian 

Motion,"  Phys .  Rev . ,  Vol.  36,  1930,  pp.  823-841. 

4.  Wang,  M.  C,  and  Unlenbeck,  G.  E.   "On  the  Theory  of  Brownian 

Motion  II,"  Rev.  Mod.  Phys. .  Vol.  17,  1945,  pp.  323-342. 

5.  Rice,  S.  0.   "Mathematical  Analysis  of  Random  Noise,"  Boll  System 

Technical  Journal,  Vol.  23,  1944,  pp.  282-332;  Vol.  24,  1945, 
pp.  46-156. 

6.  Booton,  R.  C.   "The  Analysis  of  Nonlinear  Control  Systems  with 

Random  Inputs,"  Proc .  Symposium  on  Nonlinear  Circuit  Analysis, 
Vol.  2,  1953,  pp.  369-391. 

7.  Caughey,  T.  K.   "Response  of  Van  Der  Pol's  Oscillator  to  Random 

Excitation,"  J.  Appl .  Mech . ,  Vol.  26,  1959,  pp.  345-348. 

8.  Caughey,  T.  K.   "Random  Excitation  of  a  Loaded  Nonlinear  String," 

J.  Appl.  Mech.,  Vol.  27,  1960,  pp.  575-578. 

9.  Chuang,  K. ,  and  Kazda,  L.  F.   "A  Study  of  Nonlinear  Systems  with 

Random  Inputs,"  Trans.  Am.  Inst.  Elec.  Engrs . ,  Part  II, 
Applications  and  Industry,  Vol.  78,  1959,  pp.  100-105. 

10.  Ariaratnam,  S.  T.   "Random  Vibrations  of  Nonlinear  Suspensions," 

J.  Mech.  Engrg.  Sci. ,  Vol.  2,  1960,  pp.  195-201. 

11.  Ariaratnam,  S.  T.   "Response  of  a  Loaded  Nonlinear  String  to 

Random  Excitation,"  J.  Appl.  Mech..  Vol.  29,  1962,  pp.  483-485, 

12.  Eringen,  A.  C.   "Response  of  Beams  and  Plates  to  Random  Loads," 

J.  Appl.  Mech. .  Vol.  24,  1957,  pp.  46-52. 

13.  Wiener,  N.   "Generalized  Harmonic  Analysis,"  Acta  Mathematica, 

Bd.  55,  1930,  pp.  117-258. 

14.  Ornstein,  L.  S.   "Zur  Theorie  der  Brownschen  Bewegung  fur  Systems, 

worin  mehre  Temperaturen  vorkommen,"  Zeitschrift  fur  Physik, 
Bd.  41,  1927,  pp.  848-856. 

-82- 


■83- 


15.  Van  Lear,  G.  A.,  and  Unlenbeck,  G.  E.   "Brownian  Motion  of  Strings 

and  Elastic  Rods,"  Phys.  Rev.,  Vol.  38,  1931,  pp.  1583-1598. 

16.  Press,  H. ,  and  Houboldt,  J.  C.   "Some  Applications  of  Generalized 

Harmonic  Analysis  to  Gust  Loads  on  Airplanes,"  J.  Aeronautical 
Sci..  Vol.  22,  1955,  pp.  17-26. 

17.  Thomson,  W.  T. ,  and  Barton,  M.  V.   "Response  of  Mechanical  Systems 

to  Random  Excitation,"  J.  Appl.  Mech . ,  Vol.  24,  1957,  pp.  46-52, 

18.  Lyon,  R.  H.   "Response  of  Strings  to  Random  Noise  Fields,"  J.  Acoust. 

Soc.  Am. .  Vol.  28,  1956,  pp.  391-398. 

19.  Nash,  W.  A.   "Response  of  an  Elastic  Plate  to  a  Distributed  Random 

Pressure  Characterized  by  a  Separable  Cross  Correlation," 
Tech.  Note  No.  1,  Contract  No.  AF  49(638)-328,  Engrg .  and 
Industrial  Exp.  Sta. ,  Univ.  of  Fla. ,  Gainesville,  1961. 

20.  Bogdanoff,  J.  L. ,  and  Goldberg,  J.  E.   "On  the  Bernoull i-Euler  Beam 

Theory  with  Random  Excitation,"  J.  Aero-Space  Sci. ,  Vol.  27, 
1960,  pp.  371-376. 

21.  Caughey,  T.  K.   "Response  of  a  Nonlinear  String  to  Random  Loading," 

J.  Appl.  Mech. .  Vol.  26,  1959,  pp.  341-344. 

22.  Novozhilov,  V.  V.   Theory  of  Elasticity.   English  translation,  Israel 

Program  for  Scientific  Translations,  Jerusalem,  1961. 

23.  Novozhilov,  V.  V.   Foundations  of  the  Nonlinear  Theory  of  Elasticity. 

English  translation,  Graylock  Press,  Rochester,  N,  Y. ,  1953. 

24.  Biot,  M.  A.   "Elastizitatstheorie  zweiter  Ordung  rait  Anwendungen," 

2.  a.  M.  M. ,  Bd.  20,  1940,  pp.  89-99. 

25.  Wang,  C.   "Nonlinear  Large  Deflection  Boundary  Valve  Problems  of 

Rectangular  Plates,"  N.A.C.A.   TN-1425,  1948. 

26.  Herrman,  G.   "Influence  of  Large  Amplitudes  on  Flexural  Motion  of 

Elastic  Plates,"  N.A.C.A.   TN-3578,  1956. 

27.  Tiraoshenko,  S.,  and  Woinowsky-Kreiger,  S.   Theory  of  Plates  and 

Shells.   2nd  edition,  McGraw  Hill,  New  York,  N.  Y. ,  1959. 

28.  Doob,  J.  L.   Stochastic  Processes.   John  Wiley  and  Sons,  New  York, 

N.  Y.,  1953. 


29, 


Kolmogorov,  A.  N.   Foundations  of  the  Theory  of  Probability.   2nd 
edition  of  English  translation,  Chelsea,  N.  Y.,  N.  Y. ,  1956, 


■84- 


30.  Feller,  W.   Probability  Theory  and  Its  Applications.   John  Wiley 

and  Sons,  New  York,  N.  Y. ,  1950. 

31.  Middleteon,  D.   An  Introduction  to  Statistical  Communication  Theory. 

McGraw  Hill,  New  York,  N.  Y. ,  1960. 

32.  Tolstov,  G.  P.   Fourier  Series.   English  translation,  Prentice  Hall, 

Englewood  Cliffs,  N.  J.,  1962. 

33.  Caughey,  T.  K.   "Derivation  and  Application  of  the  Fokker-Planck 

Equation  to  Discrete  Nonlinear  Dynamic  Systems  Subjected  to 
White  Random  Excitation,"  J.  Acoust .  Soc .  Am. ,  Vol.  35, 
pp.  1683-1692. 

34.  Samuels,  J.  C,  and  Eringen,  A.  C.   "Response  of  a  Simply  Supported 

Timoshenko  Beam  to  a  Purely  Random  Gaussian  Process," 
J.  Appl.  Mech.,  Vol.  25,  1958,  pp.  496-500. 

35.  Crandall,  S.  H.,  and  Yildiz,  A.   "Random  Vibrations  of  Beams," 

J.  Appl.  Mech.,  Vol.  29,  1962,  pp.  267-275. 

36.  Grobner,  W. ,  and  Hofreiter,  N.   Integraltafel .   Zweiter  Teil, 

Bestimmte  Integrale,  Springer-Verlog,  Wien,  1958. 


BIOGRAPHICAL  SKETCH 

Richard  Edgar  Herbert  was  born  November  11,  1938  at  Mount  Vernon, 
New  York.   In  June  1956  he  was  graduated  from  A.  B.  Davis  High  School. 
He  attended  the  Cooper  Union  School  of  Engineering,  a  privately 
endowed,  tuition-free  college.   Upon  receiving  the  degree  Bachelor  of 
Civil  Engineering  in  June  1960,  he  enrolled  as  a  half-time  graduate 
student  in  the  Department  of  Engineering  Mechanics  at  the  University  of 
Florida.   He  was  graduated  with  the  degree  Master  of  Science  in 
Engineering  in  February  1962.   In  September  1962  he  embarked  upon 
doctoral  graduate  work. 

The  author  was  a  half-time  Instructor  from  September  1960  to 
February  1962  in  the  Department  of  Engineering  Mechanics  except  for 
the  summer  of  1961  when  he  was  an  Assistant  in  Research.   From  February 
1962  to  September  1962  he  was  a  full-time  Research  Associate  in  the  same 
department.   When  he  started  his  doctoral  studies  he  was  granted  an 
Engineering  College  Fellowship  and  subsequently  a  Ford  Foundation 
Fellowship . 

Richard  Edgar  Herbert  is  a  member  of  Sigma  Xi  and  Omega  Delta  Phi. 


■85- 


This  dissertation  was  prepared  under  the  direction  of  the  chairman 
of  the  candidate's  supervisory  committee  and  has  been  approved  by  all 
members  of  that  committee.   It  was  submitted  to  the  Dean  of  the  College 
of  Engineering  and  to  the  Graduate  Council  and  was  approved  as  partial 
fulfillment  of  the  requirements  for  the  degree  of  Doctor  of  Philosophy. 


April  18,  1964 


SUPERVISORY  COMMITTEE 


Chairman 


■7i.£.,^.-m.^ 


— ' / — 


j/ck^tfl^^^^     cr^^-^'-t^c-^^^t_  Jt  - 


Dean,  College  of  Engineering 


Dean.,  Graduate  School 


UNIVERSITY  OF   FLORIDA 


3  1262  08553  7909 


